Resistance functions for two unequal spheres in linear flow at low Reynolds number with the Navier 

slip boundary condition 
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Resistance functions for two spherical particles with the Navier slip boundary condition in general linear 
flows, including rigid translation, rigid rotation, and strain, at low Reynolds number are derived by the method 
of reflections as well as twin multipole expansions. In the solutions, particle radii and slip lengths can be chosen 
independently. In the course of calculations, single-sphere problem with the slip boundary condition is solved 
by Lamb's general solution and the expression of multipole expansions, and Faxen's laws of force, torque, and 
stresslet for slip particle are also derived. The solutions of two-body problem are confirmed to recover the 
existing results in the no-slip limit and for the case of equal scaled slip lengths. 
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I. INTRODUCTION 

According to increasing scientific interests in micro- and 
nanofluidics and nanotechnology in recent years, fluid me- 
chanics is applied to such small-scale systems, in addition 
to molecular-level theories, where the characteristic Reynolds 
number is generally small enough to take the Stokes approx- 
imation governed by linear partial differential equations. In 
fluid mechanics, historically, both no-slip and slip boundary 
conditions were proposed in nineteenth century 1 when the 
proper boundary conditions were discussed in the first place. 
Navier 2 gave the slip boundary condition where the slip veloc- 
ity is proportional to the tangential component of the surface 
force density. For gas flows, Maxwell 3 had shown that the 
surface slip is related to the non-continuous nature of the gas 
and the slip length is proportional to the mean-free path. For 
liquids, on the other hand, from experiments at that age, the 
no-slip boundary condition was accepted and since then had 
been treated as a fundamental law. However, by recent exten- 
sive studies on the surface slip in micro and nano scales, the 
physics of the liquid-solid slip is recognized to be much more 
complicated than that for gases. Actually apparent violations 
of the no-slip boundary condition at the liquid-solid interface 
in nano scale have been reportedji 4 -^ 

Although the importance of the surface slip is realized, the- 
oretical studies and analytical solutions for the slip boundary 
condition are very limited compared with those for the no-slip 
boundary condition. Basset solved the flow of single sphere 
with slip surface, 7 Felderhof derived Faxen's law and solu- 
tions expressed by multipole expansions for single sphere 8 
and two spheres, 9 Blawzdziewicz et al. showed the interac- 
tion between the slip spheres and lubrication functions for the 
axisymmetric motion, 10 and Luo and Pozrikidis studied two 
slip spheres under the shear flow. 1 1 Recently, the present au- 
thors extended the Stokesian dynamics method (without lu- 
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brication) for slip particles using multipole expansions and 
Faxen's laws and obtained the slip dependencies for the drag 
coefficient and effective viscosity. 12 With no-slip boundary 
condition, the problem of two spherical particles is solved by 
Jeffrey and Onishi 13 and Jeffrey 14 for arbitrary size ratio of 
the particles in arbitrary linear flows. The extension to the slip 
particles was done by Ying and Peters 15 for the gas-solid sys- 
tem and by Keh and Cher. 16 for the liquid-solid system, but 
they lack the strain flows. Keh and Chen 16 applied the Navier 
slip boundary condition under a condition that the ratios of 
the slip length and radius for two particles are equal. Based on 
theory by Felderhof, 8,9 there is alternative formulation of two- 
sphere problem which covers boundary conditions of surface 
slip as well as permeability^—; they gave mobility functions 
analytically, 18 computationally, 19 and numerically, 20 and re- 
sistance functions analytically. 17 The analytical expression of 
resistance function, which is the subject of present paper, is 
limited to lower orders. 



In this paper, we will show the exact solution of two spheres 
in the form of resistance functions with arbitrary size ra- 
tio under the Navier slip boundary condition with arbitrary 
slip lengths in general linear flows including strain and shear 
flows. The present formulation is based on the no-slip case 
by Jeffrey and Onishi 13 and Jeffrey, 14 but we will show all the 
necessary equations in order that the present paper be self- 
contained. We refer equations in the references as Eq. (JO-1) 
for Jeffrey and Onishi, 13 Eq. (J-l) for Jeffrey, 14 and Eq. (KC- 
1) for Keh and Chen M 



The paper is organized as follows. In Sec|ll] the definition 
of resistance functions and Lamb's general solution are sum- 
marized. In Sec Hill the solution of single sphere with slip 
boundary condition is shown. In Sec [IV] two-body problem 
is solved by twin multipole expansions comparing with the 
results by method of reflections (shown in Appendix[A}. Con- 
cluding remarks are given in Sec [V] 
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II. FORMULAS OF THE STOKES FLOW 
A. Resistance Functions 

At low Reynolds number, the incompressible viscous fluid 
is governed by the Stokes equation 



= -Vp + pV 2 u, 
with the incompressibility condition 
V • u = 0, 



(1) 



(2) 



where p is the pressure, u is the velocity, and p is the shear 
viscosity of the fluid, Let us consider spherical particles in a 
linear flow u°° given at position x by 



u°°(x) = u°° + rr x x + E a 



x, 
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where the three constants U°°, fT, and E°° are the rigid 
translational velocity, rigid rotational velocity, and rate of 
strain of the imposed flow, respectively. According to the lin- 
earity of the Stokes equation, dynamics of the particles is com- 
pletely characterized by the resistance equation (or, equiva- 
lently, the mobility equation, that is, the inverse of the resis- 
tance equation). For two-body problem, the equation is given 
[in (J-2)] by 

U (l) - u°°(xi) 
U (2) - u°°(x 2 ) 

ra (1) -rr 
n (2) - rr 

£(2) _ £°o 

(4) 

where F s ', T (ff) , and are the force, torque, and stresslet 
of the particle a, and U^ a \ (ff) , and E (a ^ are the translational 
and angular velocities and strain of the particle a, respectively, 
and x a denotes the center of particle a. In the equation, the 
grand resistance matrix is decomposed into 6x6 submatri- 
ces. Because of the symmetry of the grand resistance matrix, 
the matrices with tilde are obtained from the counterparts as 
Bap = B^, G afi = Gp a , and H a p = W pa , (where f denotes the 
transpose) and, therefore, we need to calculate, at least, the 
rest. Following Jeffrey et al.r^^ we scale these submatrices 
[in (JO-1.7a,b,c) and (J-3a,b,c)] as 

A a p = 3n(a a + ap)A a p, (5a) 

= n{a a + a^ B a/3 , (5b) 

Cojg = n{a a + cip) C afj , (5c) 

G a p = n(a a + ap} G a p, (5d) 

H ff /j = n{a a + apf H a p, (5e) 

Mff/j = y + fl/j) 3 M ap , (5f) 

where a a is the radius of particle a, and the matrices with hat 
are dimensionless. For spherical particles, the matrices can 



be further reduced, because the geometry of the problem is 
completely characterized by the single vector r — Xp - x a . 
These submatrices are then given by scalar functions [in (JO- 
16a,b,c) and (J-4a,b,c)] as 

3f = Xipetej + Y^Sij-eiej), (6a) 

5y = Y ^ ijk e k , (6b) 

= XZw + rZiSij-etej), (6c) 

+ Y^p [eid jk + e jS/k - 2e t e je k j , (6d) 
Hijk = Y^{eiej k ie, + ej£ ik ,ei), (6e) 

M uk, = 2 ' '"t)1 ^'"tJ 

+ ~Y (ei6jie k + ej5ue k + e,5j k e\ + ejS ik ei 

-4eieje k eij 
Z M 

[SikS ji + 5 j k 5u - 6ij5 k i 

+e i e j 6 kl + 5 u e k ei + e t eje k ei 

-e,6 jie k - e jdne k - efi jk ei - e jS ik eij , (6f) 

where e = r/\r\, 6ij is Kronecker's delta, and e,j k is the 
Levi-Civita tensor. The scalar functions X, Y, and Z above 
,are called the resistance functions. We have 11 functions for 
each pair a/5. Note that for particles of other shape, such as 
spheroid for which orientation vectors should be included, the 
above factorizations by the single vector e cannot be achieved. 
From the symmetry on the exchange of particle indices a and 
P, we have the relations [in (JO-19a) - (JO-19e) and (J-5a) - 
(J-5f)] as 



*ff>M 


~ X (3-a)(_3-p)( S < ^ )' 


(7a) 




= }/ (3-ff)(3-/3)( S ' ^ )' 


(7b) 


CM 


= ~Y(3-a)(3-p)( s ^ )' 


(7c) 




= X (3-a)(_3-p)( S < ^ )' 


(7d) 




= Y(3-a)(3-p)( s > ^ )> 


(7e) 




= ~ X (3-a)(3-p)( s > A 


(7f) 




= ~Y&- a )@-p)(. s > A )> 


(7g) 


Y"p(s,X) 


= Y(3- a )(3-p)( S '^ )' 


(7h) 




= X (3-a)(_3-p)( S < ^ )' 


(7i) 




= Y(3-a)(3-p)( s > * ), 


(7j) 


^pis,A) 


~ ^0-a)(3-p)( S '^ 


(7k) 



where 



2r 



an 



(8) 



a\ + «2 a\ 

and r = \r\. Therefore, once we have obtained 22 resistance 
functions for (a/3) = (11) and (12), we can construct the grand 
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resistance matrix completely. We will see the calculations in 
Sec|IVl 



B. Lamb's General Solution 

In this paper, we utilize Lamb's general solutio n 21 ' 22 to 
solve the problem. Lamb's general solution in the exterior 
region for the pressure p and velocity u is given by 



P(r) = YjP-n-u 



(9) 



n=0 



v(r) = u(r) - M°°(r) 

oo 

= ^|Vx(^ B _,) + V<H_„_ 



n - 2 r 2 v P-n-i 



1 y f n-2 
u I 2n(2n - 



p t—i y 2n(2n - 1) p 



n + 1 p_„_i 
+ rzr 



n(2n - 1) p 



(10) 



where it 00 is the imposed velocity and v is the disturbance 
velocity field. The solid spherical harmonics p-„-u X-n-u 
and <E>_„_i are expressed [in (JO-2.3)] by 



P-n-l 



X-n- 



Pmn-{-) Y mn {6,(/>), (11a) 

m=0 

ja\ n+l 

n-l = 2^ q " m \~j Y mn{6,4>), (Hb) 



m=0 



/a\" +] 

f-n-l = Zu v mna\-J Y mn (6,(f>), (11c) 
m=0 

where F,„„ is the spherical harmonics defined by 

Y mn {e,4>) = PZ(co&e)J m *, (12) 

with the associated Legendre function P™, and p mn , q mn , and 
v mn are the coefficients to be determined from the boundary 
conditions. 



III. SINGLE SPHERE 

First, let us consider a single sphere with radius a at the 
origin. On the particle surface \r\ = a, the conventional no- 
slip boundary condition is given by 



u(r) = U + ftxr + E-r, 



(13) 



where U and J~i are the translational and rotational velocities 
of the particle, respectively. Here, we also introduce the strain 
tensor E of the particle surface, so that the boundary condi- 
tion ( fT~3b is applicable to the deformable particle at instance of 
spherical shape. For rigid spherical particle, E — 0. 



A. The Navier Slip Boundary Condition 

Navier 2 proposed the slip boundary condition, where the 
slip velocity on the surface is proportional to the tangential 
force density, as 

7 

u(r) - U + flxr + E- r + -(I- nn) ■ (cr ■ n) , (14) 

A< 

where y is the slip length, I is the unit tensor, n is the surface 
normal (equal to r/r for sphere), and cr is the stress tensor 
defined by 



cr - -pi + p [Vu + (Vit) 1 "] . 



(15) 



Rewriting Eq. ( fT~4T > by using the disturbance field v and the 
imposed flow u°°, we have 



V 1 = w + —t 



(16) 



where the disturbance part t and imposed part t°° of the tan- 
gential force density are defined by 

t = (I — nn) ■ (cr Y ■ n) , (17a) 
t°° = (I - nn) ■ (cr°° ■ n) , (17b) 

and the corresponding stresses are 

cr v = -pi + fi [Vv + ( Vi)) f ] , (18a) 
<t°° = p [Vu°° + (Vm 00 ) 1 "] , (18b) 

w A is defined by 

w A = AU + Aft x r + AE ■ r, (19) 

and AU = U - U°°, Afl = fl - fl°°, and AE - E - E™. 
From the imposed flow in Eq. (O, t°° becomes 

t°° = ^ (J - nn) ■ E 00 ■ r. (20) 
r 

Note that, in the slip boundary condition ( [Tol l, the left-hand 
side is the disturbance quantities and the right-hand side is the 
imposed quantities. Also note that, on the imposed part, the 
slip contribution appears only on the flow with E°° + as 
shown in Eq. (f2Qb . 

In terms of Lamb's general solution for the disturbance field 
v in Eq. ([Tol l, the corresponding surface force density / is 
given b y 21 ' 22 

/ = cr v -n 

= ^J]l-(n + 2)Vx(rx-n-i) 

n ^ 

-2(n + 2)V0)_„_ 1 
l(n+l)(n-l) 9 _ 
p n(2n - 1) 

1 2« 2 + 1 

7^ ttT-«-i 

p n(Zn — 1) 



(21) 
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and t defined in Eq. ( I17al i is expressed by 

* = - r Y\- {n+2)Vx{rx -"-' ) 



-2(n + 2)|v-~W_„_ 1 

l (n + lXn-l) 2 / r3 , 

+ 7^ tt— r V - -— /?_„_] j>. (22) 

n(2n - 1) \ r or 



1. Three Scalar Functions 

In order to achieve the boundary condition for Lamb's gen- 
eral solutions, Jeffrey and Onishi 13 used three scalar functions 
as in Happel and Brenner, 22 §3.2. Consider a general vector 
field g and its surface vectors G defined by 



Because v satisfies V ■ v = 0, Vdiv is given by 



Vdiv = -rV ■ V = r—v r \ , 

Or \\r\=a 



(29b) 



from Eq. ( [27) . V ra d is independent of its radial component V r 
as shown in Eq. ( 126bb , so that it is simply written by v as 



Vmt = r ■ V x V = r ■ V x v\ 



l|r|-a 



(29c) 



From Lamb's general solution for v in Eq. ( fTOb . then, the 
three scalars are obtained as in Jeffrey and Onishi. 13 

b. Tangential Surface Force Next, let us consider t 
which is necessary for the slip boundary condition ([T6V Its 
surface vector is defined by 



G(0,<f>)=g\ 



\r\=a 



so that 



dG 
~dr 



= 0. 



We define the following three scalar functions 



G ra d 
Gdiv 



- • G, 
r 

-rV G, 
r-VxG. 



(23) 



(24) 



(25a) 

(25b) 
(25c) 



T(0, 4>) = t\ 



r\-a 



The radial component of t is zero by definition as 



- T = 0. 



From Eq. d27] i. therefore, we have 



r div = -rV • T = -rV ■ t\ 



\\r\=a 



(30) 



(31a) 



(31b) 



Obviously, the first scalar G m d is the radial component G r = 
(r/r)-G itself. The other two, Gd; v and G rot , are related to the 
tangential components (i.e., Gg and in polar coordinates), 
except for the factor -2G r on the divergence, as 



Gdiv 

G rad = y 



1 06 sin6»/ sin6» d<p 



J_dGg + ld_ + cosfl\ 
"^me~dS + \b~6 + sine 



(26a) 
(26b) 



where is +1 in the right-handed coordinates and -1 in the 
left-handed coordinates. It should be noted that the divergence 
of the surface vector G is related to the 3D vector field g as 



I d I 

Gdw =-rv-g\ +r— g r , 

\\r\=a Or \\r\=a 



(27) 



where the substitution of \r\ — a is applied after the deriva- 
tives. 

a. Velocity Field As a first example, consider the distur- 
bance velocity v, whose surface vector is defined by V as 



V(G, </>) = v\ 



\\r\=a 



By definition, the first scalar V m d is given by v as 



r r i 

= - ■ V = - ■ v . 

r r i\r\=a 



(28) 



(29a) 



Because the rotation has no radial component for an arbitrary 
vector field, we can use the bare surface force / for the bound- 
ary condition for the tangential force t as 



T mt =r -V xT = r -V xt = r • V x / . (31c) 

l|r|=<7 l|r|=fl 



Using Lamb's general solution in Eq. (l22l i. the three scalar 
components for t are given by 



-tt = 0, 



-rV-f - -a* J) 



2n(n + l)(n + 2) 



(n + l) 2 (n - 1) 



2n- 1 



(32a) 



(32b) 



r-Vxt = -^J](ra + 2Mra + l)^_„_i. (32c) 



c. Disturbance Part Three scalars for V are obtained 
by Eqs. d29a| ), ( |29b) , and d29c| ), and the slip contribution 
-(ylfi)T by Eqs. d32al >. d32bl . and d32cb . Substituting Lamb's 
solution ( [Tol l with the expansions in Eqs. (II lab , (II lbl ). and 
(II lcb and putting r — a, the three scalars of the disturbance 
part, i.e. the left-hand side, of the slip boundary condition 
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([TBI are given by 



2. Recurrence Relations 



v- 7 -t) 



oo n 



n + 1 



/ div 



oo « 

EE 

w=0 m=0 



2(2m-1) j 
(«+ l)(n + 2)(l +2ny)v m 



(33a) 



«(n+ 1) / 2(«+ l)(n- 1), 

1 + 7\Pr, 



\v- y -T\ 



2(2n-l) 

Fm„(0,0), (33b) 

CO /J 

\ #~ / rot n= o m =0 

(33c) 

where the scaled slip length ■y is defined by 

y 

a 

of. Imposed Part Let us look at the three components for 
the vector u> A in Eq. ( fT9] l. Note that the divergence is zero as 
shown by 

diwf = ^tAfljtf a + AEytfy = 0, (35) 

because E& = 0. Therefore, we need to calculate the diver- 
gence component through the derivative of the radial velocity 
(as for v). The three components for w A are then given by 



r 



(34) 



r '„,A 



n 



-wf = -AUi + —AEa, 
r r r 

rieijkdjW^ = 2r ! Afi i . 



(36a) 

(36b) 
(36c) 



We use the identity £ijk£jki = 25 ,7 for the last equation. For 
the three components are given as follows. The normal 
component is zero by definition as 



-<T = 0, 

r 



(37) 



and, therefore, the divergence component is obtained through 
Eq. ( |27| | as 

- rdfi = -2/irdi [dtjj - ^jf-)Ej k = 6^EJ k , (38) 
where we use E^ k = 0. The rotation vanishes as 

mjkdjtt = tymjkdj [sj-j - r -^^E2 = 0. (39) 

Define the surface vector of the right-hand side of the slip 
boundary condition ( [TBI by 



w = U A + -*° 



\r\=a 

The three scalars for W are then given by 
W m d = e,A{7j + ejejaAEjj, 
Wdiv = eiejaAEij + GyeiejaE^j, 
W rot = 2e / aAQ / . 



(40) 



(41a) 
(41b) 
(41c) 



Let us introduce the spherical harmonics expansion for the 
three components of the imposed part by 

oo n 

w rdd = YjTj XmnYmn(e '®' (42a) 

«=() m=0 
00 n 

Wd > v = YjTj*™ 7 ™^®' (42b) 

11=0 m=0 

00 n 

W mt = YjYu *™*™^®- (42c) 

«-0 m-0 

From Eqs. (141 al l. (141bb . and (141cl >. the coefficients ^ m „, ^ m „, 
and oi m „ are given by the parameters AU, Af2, AE, and 
£7°°. Therefore, by the boundary condition ([TBI at the surface 
|r| = a with the scalars of the disturbance fields in Eqs. (133al i. 
( 133bl ). and ( 133cl ). the coefficients (p„ m ,qmn,v mn ) are given by 
the boundary condition (^ m „, <^ m „, w m „) as 

2n ~ 1 r 



n + 1 
(« + 2)(2« - 1) 



n + 1 



^2n,2n+lXmn, (43 a) 



1 



2(71 + 1) 

n 



(43b) 
(43c) 



2(n+ 1) 
1 

1 mn ~ „/.. , ^ ^0,114-2^, 



where 



n(n + 1) 



r 

1 m,« 



1 + my 
I + try 



(44) 



Note that in the no-slip (y — 0) and perfect-slip (y = oo) 
limits, r„, „ reduces to 



r - 



|1 for y = 0, 
Im/n for y = oo. 



(45) 



B. Single Body Solutions 

In the following, we solve single-body problem with the 
slip boundary condition through Eqs. (143a) . (143bl i. and (I43cb . 



1. Translating Sphere 

Consider translating sphere with the velocity U - (0, 0, U), 
which is given by 



Xm.n — U6o m 5i n . 



(46) 
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Substituting the condition d46b into the recurrence relations 
(|43a| ), d43b| >, and d43c| ), we have the solution 



Pmn - -jUT2,z5 m o5 n \, 

v mn = ^LT 0; 3(5 m 0(5„i, 
tfmn — 0. 



(47a) 

(47b) 
(47c) 



The force acting on the particle is given by the coefficients of 
Lamb's general solution [in (JO-2.10)] as 



F = \npa [p m z - p u (x + iy)] . 



(48) 



where x, y, and z are the unit vectors in x, y, and z directions, 
respectively. Therefore, the force on the sphere translating 
with the velocity U in z direction is 



F = 6nfxaT2,3Uz. 



(49) 



This is identical to the result by Basset. 7 (See also Lamb^i Art. 
337, 3° and Felderhof.) 8 Substituting the coefficients J47al >. 
d47b| ), and d47c| i into Lamb's general solution in Eq. ( TTOb 
and rewriting the parameter U by the strength of the force 
F through Eq. ( |49l , the disturbance field is given by 

^i( 1+r «4 v2 ) J - F ' (50) 

where J is the Oseen-Burgers tensor 

W=\{6ij+ r -^\ (51) 

2. Rotating Sphere 

For the problem of rotating sphere, W mt in Eq. ( 14 let is the 
only non-zero component. Consider a sphere with the angular 
velocity Q - (0, 0, Q), which reduces to 



2aQ,So m 5i„- 



(52) 



Substituting the condition d46b into the recurrence relations 
( I43al ). ( 143bl >. and d43cl ). we have the solution 



Pmn = 0, 
Vmn = 0, 

q mn = aQTo^6 m o6„i. 



(53a) 
(53b) 
(53c) 



The torque acting on the particle is given by the coefficients 
of Lamb's general solution [in (JO-2.11)] as 



T = Snpa 2 [q m z -q u (x + iy)} . 



(54) 



Therefore, the torque on the sphere rotating with the angular 
velocity Q in z direction is 



T = 87rpa 3 r ll3 Clz. 



(55) 



This is consistent with the result by Felderhof 8 and Padma- 
vathi et al— Note that the torque T would vanish for the 



sphere with the perfect-slip surface (for y = oo). Substitut- 
ing coefficients d53at , d53b| ), and d53c| ) into Lamb's general 
solution in Eq. ( TTOb and using Eq. d55t . the disturbance field 
is given by 



1 



-R T, 



where 



Rijir) = e ijk -j. 



3. Sphere in Strain Flow 



(56) 



(57) 



For the problem of sphere in strain flow, we have two non- 
zero components on W. Here we assume the rigid sphere, so 
that E = and from Eqs. ( 141 at , d41bt , and d41c| >. 



W md = -eiejaE™-, (58a) 
W div = -e,ejaE™(l-6y), (58b) 
W mt = 0. (58c) 



Let us consider the strain given by 



This is achieved by 



(59) 

(60a) 
(60b) 



Substituting the boundary conditions d60ab and (I60bb into the 
recurrence relations (I43at . (143b) . and (!43cK we have the solu- 
tion 



Xm,n - -aESo m O~2n, 
2 

i/r,„,„ = -aE (I - 6f) 6 0m 6 2n . 



Pmn - —aET2,5d m Q5 n 2, 

v mn = -^ a Er 5S m o6„2, 
qmn — 0. 



(61a) 

(61b) 
(61c) 



The stresslet acting on the particle is given by the coefficients 
of Lamb's general solution [in (J-6)] as 

S = 2npa 2 |/?02 {zz - ^7 

-p n [xz + zx + i (yz + zy)] 
+2p 22 [xx-yy + i(xy + yx)]}. (62) 

Therefore, the stresslet on the sphere in the strain flow with 
the parameter E is 



S = ^7rpa 3 r 2 , 5 E, 



(63) 



which is identical to the result by Felderhof.- This yields to 
the effective viscosity jj* of the suspension in the dilute limit 
up to 0((f>) as 



u* 5 

- = i + -r 2 , 5 « 

Lt 2 



(64) 



where <p is the volume fraction. This is identical to the expres- 
sion (9-5.11) in Happel and Brenner. 22 The effective viscosity 
of slip particles has two extremes as 



yu* 1 1 + j(p for no-slip particles, 
H (1+0 for perfect-slip particles. 



(65) 



The latter agrees with the result for spherical gas bubbles. 
Substituting the coefficients ( 161 at , d61bt , and ( 16 let into 
Lamb's general solution in Eq. (TlOb and using Eq. d63l >. the 
disturbance field is given by 



where 



K ijk {r) = -3 



(66) 



(67) 



IV. TWO-BODY PROBLEM 

Now, we study two-body problem. We will determine 22 
resistance functions mentioned in Sec. Ill Al Following Jeffrey 
et al.r^^ we write these functions in terms of the coefficients 
f m and determine the coefficients. Here we summarize the 
definitions of the coefficients: X^„ are given [in (JO-3.13) and 
(JO-3.14)]by 



X' 



m fXA 



x 



m=0,even 

-2 



■XA 



m=l,odd 



(68a) 
(68b) 



Y* p [in (JO-4.13) and (JO-4.14)] by 



m=0,even 

-2 



J m 



E 



m=l,odd 



(69a) 
(69b) 



[in (JO-5.3) and (JO-5.4)] by 
Y'Ms.A) = J] 



[(l + A)*]"' 



yf 2 (M) 



~* y Jm 

(1 + A)2 Zj [(1+^if 

m=0,even 



(70a) 
(70b) 



X^ [in (JO-6.7) and (JO-6.8)] by 



X 



m=0,even 



AT 



2 = (TT? f dd [0T 

m=l,odd 

F ( C [in (JO-7.7) and (JO-7.8)] by 

^(M) = X ^ 



[(l+W 



o °° fYC 

V C I X\ - V -An 



X° [in(J-18a,b)] by 



X 



.1 111 



^ ' Zj [(1 + 



[(1+A)s]" 



vg m = — 1_ y _^ 

12 (i +^) 2 v [(i +. 

m=2,even 

7« [in (J-26a,b)] by 

fYC 

.) Ill 



iii(M) = J] 



m= 1 ,odd 

-4 



Ki + ii)*]- 



m=0,even 



7" [in (J-34a,b)] by 



rFfl 



«m - E 

m=0,even 



X* [in (J-47a,b)] by 



X 



m=u,even 



x 



7$ [in (J-63a,b)] by 



[(i +A) s ] m ' 



*<•■* - E nfw- 

m=u,even 



f 

.1 in 



8 



[in (J-78a,b)] by 



r l = r 3-a + r 2 - 2r 3 _ Q .rcos6'3- Q ,. 



(83b) 



f z 

■' 111 



Z 



fZM 



m=\ ,odd 



(78a) 
(78b) 



A. Twin Multipole Expansions 

Let us consider two particles a = 1 and 2, whose centers, 
radii, and slip lengths are given by x a , a a , and y a , respec- 
tively. The scaled slip length for particle a is defined by 



- Ja 

y a = —■ 

Q,(y 



(79) 



First, we outline the derivation of equations among coeffi- 
cients (p mn , q m „, v m „) and (iffmn,Xmn, w») for the slip spheres. 
Then, we solve the recurrence relations for each problem and 
obtain all the resistance functions. 



1. Outline 

In Sec. |IIIJ the problem of single slip sphere has been 
solved by Lamb's general solution ( TTOb through three scalars 
of the surface vector on the both sides of the slip boundary 
condition ( TToT l. Jeffrey et alM^ solved two-sphere problem 
with no-slip boundary condition, i.e. y — in Eq. ( TToT l. To 
complete the boundary condition for two slip spheres, we need 
to obtain the tangential force density caused by particle (3 - a) 



on the surface of particle a. Let us denote it by i' , that is 



t' ia) = (l - n ( V) • ( 



P (3-oD 



(80) 



where n (t?) is the surface normal of particle a (r^ a) /r a for a 
sphere), and er (3 ~ Q,) is the disturbance part of the stress caused 
by particle (3 - a) given by 



r (3-a) 



,<?-*> j + 



(81) 



Here, p (3_Q,) and v®~ a ' are expressed in terms of Lamb's gen- 
eral solution for the polar coordinates of particle (3 - a) given 
by Eqs. © and ([Tol l, respectively. Because er (3 ~ ff) • n (a) + 
/ (3 ~ tv) , we cannot use the surface force density / in Eq. ( |2TI ). 

Following similar calculations by Jeffrey and Onishi 13 for 
the disturbance velocity, we can write <r (3 ~ ff) by the spherical 
harmonics with respect to the particle a in terms of the trans- 
formation [in (JO-2.1)] 



' s=m * ' 



and the following relations [in (JO-2.7)] 



V a,4>) 

(82) 



TV = 



h-a (ri- a - f cos 6» 3 _„) + 3 _ Q r sin6» 3 _ Q 



(83a) 



After substituting the expansions for the solid spherical har- 
monics PL~"l, X--i> an d $>_"l m Eqs. ( II lab . ( II 1M . and 
( II let , the three scalars of the surface vector of t'^ a ' are ob- 
tained in the form of the expansion with spherical harmonics 
Y mn (8 a ,(p). Combining the results of t' (a) and those of the 
disturbance velocity on particle a caused by particle (3 - a) 
given by Jeffrey and Onishi 13 with the single-sphere problem 
in Eqs. ( I43al ). ( I43bl ). and ( I43cb . we have three equations for 
the coefficients, corresponding to Eqs. (JO-2.9a), (JO-2.9b), 
and (JO-2.9c) for the no-slip case, as 



^-(n-Y){\-2(n+Y)y a ) X \ 



(n + l)(2n+l)(l+2r tt )v; 



(a) 
inn 



n+l 



p (B) 

i inn 



(84a) 



^+(n+2)(l+2n%)x 
n+l 



Ur) 
inn 



2/7 



^ r (i + (2« + Dr ff )ps + 2(; + + Mc 1 ^ 

s=m \ ' 

x [i(-l) a m(2n + 1) (1 + Zy a )q%; a) t 3 - a 

+n(2n +1)(1+ 2y a ) v^ a h\_ a 

2n + 1 , _ . 

+ 2n^l {1+2 ^ 

ns(n + s — 2ns — 2) — m 2 (2ns — 4s - 4n + 2) n_ a ^ 
x ~ ~ — ; p ms 



2s(2s -!)(« + s) 



n 

-2 1 

» 



,(3-«) f 2 



(84b) 



Jo) 
i inn 



n{n + l){\+{n + 2)y a )ql 

s=m * ' 



(3-o) 
ms 



where 



- Ua 
t a — 

r 



, (84c) 



(85) 



In writing these three equations, we can take any independent 
linear combinations in principle. For the single-sphere prob- 
lem, we may write three equations for p mn , q mn , and q mn as in 
Eqs. ( I43al i. ( I43bl i. and ( I43cl i. or those for the coefficients of the 
boundary condition Xmn, ftmn> an d w,„„, instead. Jeffrey and 
Onishil 3 - take equations for i/>^ - (n - l)x ( mL ifftH + {n + 2)x ( ml 
and a>„ a ^ for no-slip particles. Here we extend the equations 
for slip particles so that the structures of the equations for no- 
slip case would hold, that is, the interaction terms (with the 

, (3-af) . 



summation of s) contain only p mn in Eq. ( I84al ). and the term 
of is eliminated in Eq. J84b| ). Equation d84c| ) for <J^ n is 
just the same choice to the no-slip case. 
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Note that Keh and Cher. 16 take a different form for the first 
equation, that is, ifr^-iXn - 1) + (2n 2 + l)y a )xmn in Eq. (KC- 
20a). Although they are mathematically equivalent, Eq. d84at 
is simpler and we will use it later in this paper. Also note that 
there are typos in Keh and Chen 16 at Eqs. (KC-20a,b,c) where 
P(3-a) (T(3-«) m P resent notations) should be replaced by f}~A, 
If we look at the slip boundary condition from which these 
three equations are derived, it is obvious that only the slip 
length of particle a would appear there. It should be noted 
that the results such as coefficients fa in Keh and Client are 
correct, because they took a simplification that the scaled slip 
lengths for two particles are the same (yi = 72 m present 
notations). 



2. Recurrence Relations 

For resistance functions, the boundary conditions are given 
completely by Xmn> tymn, an d u>„ m , which are independent of 
the distance between the particle r and therefore t a and t^ a . 
This means that the coefficients P npq , V npq , and Q„ pq of the 
(p, g)-term in the expansion by ^t\_ (see, for example, Eqs. 
d89ab and d89bb in the following) are solved by the recurrence 
relations for p > and q > with the initial condition at p = 
and q — 0. Therefore, we split the above three equations into 
two parts, the initial conditions and the recurrence relations. 
The initial conditions are 



Imn 



= Z 



-i(-D c 



(B+l) 



(«-l),n+2^ 



(3 -a) J! 



ar'3-ar 



ns(n + 1) 



-(n-l),n+ 



O-a) t n f s 



(87c) 



It should be noted that the initial conditions are independent 
of m, while the recurrence relations are not. Therefore, the 
initial conditions are the same for X (m = 0), Y (m =1), and 
Z (m = 2) functions for each problems (translating, rotating, 
or in the strain flow). 

We also note that the recurrence relations have ar-dependent 
quantity r (t?) , so that we need to solve the coefficients P npq , 
V npq , and Q npq for a as well as (3 - a), while, for the no-slip 
case, the coefficients for a and (3 - a) are identical. 

The results shown in the following are obtained by the pro- 
gram implemented on an open source computer algebra sys- 
tem called "Maxima". 24 The program is relatively slow due to 
its symbolic calculation and the coefficients are obtained up to 
k = 20, at least. We also implement a code in C with floating- 
point variables where the parameters a a and y a must be given 
by numbers for the calculation. With this code, we can obtain 
the coefficients around k = 100. 



1 nm 



In- 1 



n + 1 

(n + 2)(2n 



1 0,2«+lV»m 



2(2n + l)v 



(n) 
inn 



+ !2r w v 

H + 1 2n,2n+l A -- 

n + 1 °- 2Wmn (n+1) 



(a) 
inn ' 



T<a) 



(86a) 

(a) 



I run 



+r (a) D (a) 

0,2 "mn ' 
1 



» 



n(n + l) °'" +2 
The recurrence relations are 



(86b) 
(86c) 



I mn 



Zl n + s \ 
\n + m) 



-i(-l) ff m 



(2ra+ l)(2n- 1) 
n+1 



*~ 1 I — 1 1 *4 m p * * 1 _ ™ 



- 2,2«+l l !m.s 



,.(3-a)/i-l f *+2 
2,2n+l ms l a ( 3-a 



n(2n + l)(2n - 1) 
n + 1 

2n + 1 ns(n + s — 2ns - 2) — m 2 (2ns — 4s — 4n + 2) 
n + 1 



2s(2s - l)(n + s) 



<<*) _(3-a)j»-l f j 
2,2n+\"ms l a lr h-a 



n(2n- 1) (a) 
' 2.n - I ) "•"• | 



(3-a)jj+l„i 

f fill! 



3-a 



(87a) 



B. X Functions (m = 0) 

For the case of m = 0, and g <3 ~ 1!,) are decoupled from 
the others. 



i. X 4 Function 



The boundary condition for the X problem is given by 



X%> = U6 m0 S n u <C = 0, < J =0. (88) 

To obtain the coefficients for each order of the power of r, we 
expand the coefficients [in (JO-3.4) and (JO-3.5)] as 



(a) 



CO OO 



P =0 q=0 

CO CO 

zz 



(a) 
npq 



2 ^^ 2 ( 2 " + D 



(89a) 
(89b) 



2(2n + IV 



r(ff) n W 

1 0,2' mn 



-z 



2n 



(87b) 



n + m j ( n + \)(2n + 3) ~ (2 " +1) 



(3-ory+l. 



3-q" 



Substituting the expansions, we have the initial conditions for 
p — and q — from Eqs. ( I86ab and (I86bb by 



p( a ) _ r w '\A a> — a r w 

^«00 ~ ""l 1 2,3' v n00 _ 0,111 0,3' 



<«) x/or) 



(90) 
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and the recurrence relations for p > and q > from Eqs. 
( I87al > and ( fSTbb by 



p(a) 
«P9 



z 



n + s 
n 



n(2n - l)(2n + 1) 



.(<!) 



V 



(3-a) 



2(n + l)(2s + 1) 2 - 2n+1 «(?-*-2)Cp-»+D 



n(2n + 1)(« + 5 — 2ni — 2) 
2(n+ l)(n+ s)(2s- 1) 



3 (3-o) 



«(2n - l)_ r („) 



2{n + 1) 



0,2«+l ,s(g-s)(p-n-l) 



2,2n+l s(#-s)(p-n+l) 

(91a) 



y(a) 
r npq 



-z 



ria) p(a) 
1 0,2 r npq 



n + S 



(91b) 



2n 



(n+ l)(2n + 3) -( 2 " +1 >- 2 



3 (3-Q-) 
s(«-s)(p-n-l)' 



The initial conditions correspond to Eqs. (KC-26a,b) and the 
recurrence relations to Eqs. (KC-27a,b). Note that Eq. ( 191bt 



for Vnpq is simpler than the corresponding equation in Keh 
and Chan (KC-27b), because we use the simpler recurrence 
relation in Eq. (|84a| >. 



The coefficient tf Aa is defined [in (JO-3.15)] as 



pends on a. The explicit forms up to k = 1 are 



fS 

ff 



■XA\ 



pXA\ 



cXAl 



cXA\ 



ff 



- (iS). 
= ^K'rg), 

= ^(9 ( r«) 2 rg), 

= 4-<3 r S) 
+ A 2 (27(r«) 2 (r®) 2 ) 

+ ^Higig). 

= A(-24lglglg) 

+ A 2 (8i(r^) 3 ( rg) 2 ) 

+ A 3 (l2(T«) 2 (5rg-2r®)), 

= A 2 (36r« ( rg) 2 (5r^-3r^)) 
+ ^(243 ( r«) 3 ( rg) 3 ) 

+ * (36ffS)>lS(5lg-3lg)). 

= ^(i6 ( r<:») 2 rg) 

+ ^(i08(4;>) 2 (rg) 2 ( 5 r^-<>)) 



+ A 



(-r^(48or; 



_ 32r (D r (2) ) \ 
JZjL 0.3 1 0,3 >) 



<l) r (2) 
0,3 2,5 



729(rg) 3 (rg) 3 



■(2) 
3 



(93a) 
(93b) 
(93c) 



(93d) 



(93e) 



(93f) 



(93g) 



+ A 4 (2i6(r«) 3 rg(5rg-2rg)) 
+ ^ 5 (i6(r^) 2 (i26r^ - 9org + 5r®r| 

44I^ 2 )/5), 

= ^ 2 (48(rg) 2 ( i26rWr^ - 70i$ig 

- 45r o! 5 r 2,3 + 15 2 + 4r ? 3 j)/5) 
+ ^ 3 (i620(r^) 2 ( rg) 3 (2r^-r^)) 

+ A " ( 3r 2!3 r S( 800r 2 ! 5 r 2 2 5 " 560r o!3 r 2 2 5 

-56orgrg + 729(rg) 3 (rg) 3 + 96r$r£>)) 
+ a 5 (i620(r^) 3 (rg) 2 (2rg - r®)) 
+ a 6 (48(r^) 2 (i26rgrg - 7or®r® - 45rgrg 

+ 15(rg) 2 + 4r (2 ( 3 )/5). (93h) 



^ = 2^%^. (92) 

9=0 



Here we see a slight difference from the no-slip case. This is 
because of the a dependence of P " pq , so that f* Aa also de- 



The results are identical to those obtained by method of re- 
flections in Eqs. (IA14cb . ( IA14al ). and dA14bl ) for the terms 
containing one or two T's, because only the first reflection 
from the particles 1 to 2 is taken and the higher reflections 
are missing in the present calculation of the method of reflec- 
tions. Therefore, f* Al and A 2 term in f^ A do not appear. 

Also the results reduce to those by Jeffrey and Onishii 3 - in 
the no-slip limit 7 = 0, and those by Keh and Chen 16 in the 
case of yi = 72. Therefore, they also reduce to those by Het- 
sroni and Haber 25 in the perfect slip limit 7 = 00. 
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2. X a Function 

With the same recurrence relations and the initial condition 
for X A , that is, for the translating particles, the function X G is 
obtained from the coefficient Pi pq for the stresslet instead of 
P\pq for the force. 

In this case, the coefficient f* Ga is defined as 



(94) 



The explicit forms up to k — 7 are 

XCl _ n 



fi 





fXGl 



(95a) 
(95b) 
(95c) 
(95d) 



f r = ^(isrgrg), 

/f* = A(45I<SlglS), 

jf* = A(-3fflgrg) 

+ A 2 (l35r«(^) 2 4») 

+ A 3 (-60rgrg), (95e) 

/f* = A(-i2rg(5r«rg + 9 r«4»)) 
+ ^ 2 (405(r£) 2 ( rg) 2 r£) 

+ ^ 3 (i2or^r^(5rg-2rg)), (95f) 

= ^(36(^)2(25(1^)*- 101^1^-91^1^)) 

+ A 3 (i2i5(rg) 2 (rg) 3 rm) 
+ ^ 4 (90or^rgr«(2rg-r®)), (95 g ) 
jf» = A(i44r£>r«rf>) 

+ ^ (losr^rg) 2 ^') 2 - i 5 r<»r« - 9igrg)) 
+ ^ 3 (-3(80or$rgrf> + 96or«r«F 



-(2) 
2.5 



-I2i5(4») 3 ( rg) 3 r« 

" 8Ui 0.3 1 0.3 1 2,5 



■48rgr«r«)) 
+ /(i620 ( r^) 2 rgr«(5rg-2rg)) 
+ A 5 (48r«r«(i26rg-9org + 5 rgrg 

+4r%)). (95h) 

The results are identical to those obtained by method of re- 
flections in Eqs. ( IA20al ). ( IA20bl ). and ( IA20cl ) for the terms 
containing one or two F's, similarly to X A . The results reduce 
to those by Jeffrey 14 in the no-slip limit y = 0. 

3. X c Function 



Note that Q npq is decoupled from P npq and V npq for m = 0. 
Using the expansion [in (JO-6.4)] 



y o« Zj Zj y »w " 

p=0 ?=o 



oo oo 



(97) 



we have the initial condition for p — and g = from Eq. 
(l86cl by 



QnOO - 5„ir| 



(«) 
0.3' 



(98) 



and the recurrence relation for p > and g > from Eq. J87cl > 
by 

0(<U = V ( " + 5 \ _L_ r w (3 - ff) (99) 



The coefficient is defined as 



(100) 



9=0 



where j - for even and y = 1 for odd k. Because many 
terms of f* cl in lower orders are zero, we show the explicit 
forms up to k — 1 1 as 



f xc\ 
Jo 




K). 


(101a) 


fXCl 

J\ 




0, 


(101b) 


fXCl 
J2 




0. 


(101c) 


fXCl 

J3 






(lOld) 


fXCl 

■h 




o, 


(lOle) 


fXCl 
J5 




0, 


(lOlf) 


f XCl 
J6 




A 3 (64 ( rg)2rg), 


(101g) 


fXCl 

■h 




o, 


(lOlh) 


fXCl 

J% 




^(768r (2 1 ) 4 (r^) 2 ), 


(1011) 


fXCl 
J9 




/(5i2 ( r^) 2 ( rg) 2 ), 


(101j) 


fXCl 

J 10 






(101k) 


fXCl 

Ju 




/(6i44r^ 4 r^ ( rg) 2 ) 






+ 


A 8 (6i44^ 4 (rg) 2 rg). 


(1011) 



The results are identical to those obtained by method of re- 
flections in Eqs. ( lA44al i. ( IA44bl ). and ( IA44cl ) for the terms 
containing one or two F's. The results reduce to those by Jef- 
frey and Onishi 13 in the no-slip limit y - and those by Keh 
and Chen 16 in the case of yi = 72- 



The function X c gives the torque for the rotating particles 
in the axisymmetric case (to = 0). The boundary condition is 
given by 



y (a> 



= o, <^ = o, 



(0, 



- 2U5 m ()6 n \ 



(96) 



4. X M Function 

The function X M gives the stresslet under a strain flow in 
the axisymmetric case (to = 0). Therefore, it is derived by the 
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coefficient P% pq for the stresslet from the same recurrence re- 
lations for X A with a different initial condition. The boundary 
condition is given by 



C. Y Functions (m = 1) 

1. Y A Functions 



Xmn = -^a a E a 6 Qm 52n, (102a) 
0jg = ~a a E a (l-6y)d 0m d 2 n, (102b) 

wSS = 0' ("He) 

which corresponds to Eq. (J-41) with the correction due to the 
slip. Using the expansion [in (J-42) and (J-43)] 



Pan ~ a a aEg ^ i ^ i Pj,}„t„t. 



10 

T 



oo oo 



(103a) 



p=0 g=0 



(a) 
"(to 



V 



(a) 



/7=0 9=0 



^— iS/f _. (103b) 



the initial conditions for P npq and V npq are given from Eqs. 
(I86at and < [86bl by 



0,5 - 



(104) 



The coefficient A is defined as 



The boundary condition for the Y A problem is given by 
X™ = {-\) a US m ,8 n i, ^ = 0, w£» = 0. (107) 



mn 

,i 13:,. 



(Note that the equation by Jeffrey and Onishi, in p. 271, lost 
the factor U for;^.) Again, we expand the coefficients by t p a 
and d as 

3-a 



(a) 
Pin 



(a) 
Tin 



CO CO 



p=0 q=0 
oo oo T,(a) 

( L) 2^ZjZj2(2» + 1) a 3 -' 

p=0 g=0 v 7 



CO CO 



p=0 9=0 



(108a) 
(108b) 
(108c) 



Also note that the minus sign in the right-hand side of (JO- 
4.5) is missing. Substituting these expansions into Eqs. ( I86al ). 
( I86bl ). and ( I86ct . the initial conditions are given by 



fXMa _ r,k V" 1 p(a) 
9=0 



(105) 



where j — for even k and /' = 1 for odd k. The explicit forms 
up to k — 7 are 



/<r = (iS). 



(106a) 
(106b) 
(106c) 
(106d) 
(106e) 



/ 2 AM1 = 0, 

/P" = ^ 3 (4or^rg), 

/r 1 = ^(6org ( rg) 2 ), 

/™ = ^(-192^^) 

+ ^(isorgrgrgrg) 

+ ^ 5 (-i92r®r^), (i06f) 

jf™ = A(-288lglglg) 

+ ^(54org ( rg ) 2 (rg) 2 ) 

+ ^(i60(r^) 2 (iorg-3r®)), (i06g) 

/f" = / (48rg(50(r^) 2 - 2or$r™ - 9r^r^)rg) 

+ ^ 5 (i620(r^) 2 (rg) 2 r^rg) 

+ /i 6 (48r^r^(50(rg) 2 - 2orgrg 

- 9r o 2 5 r 23))- ( 106h ) 



The results are identical to those obtained by method of reflec- 
tions in Eqs. ( IA62ab , (IA62bt . and (IA62ct for the terms con- 
taining one or two F's. The results reduce to those by Jeffrey 14 
in the no-slip limit y = 0. 
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nl 1 2,3' 
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nOO 



.(«) 

0,3' 



o ia) 



0. 



(109) 



which correspond to (KC-37a,b,c). From Eqs. ( I87at . (187bl >. 
and (I87cb . the recurrence relations are given by 



npq 



§(::;) 
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n + 1 
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»(2»- l) rW (3 - ff ) 
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y(a) _ p(«)p(ff) 
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(110a) 
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y / n + s\ In (g) (3 _ ff) 
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Note that Eqs. (Ill Oat and (111 Ocb correspond to (KC-38a) and 
(KC-38b), while Eq. ( fTTObb is simpler than Eq. (KC-38c). 
The coefficient fY Aa is defined as 



cYAo 



(111) 



q=0 



The explicit forms up to k — 7 are 
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- (iS). 

= *(3l£lg/2), 

+ i 2 (27(r^) 2 (rg) 2 /s) 
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(2) r (l)\ 
0,3 2,3/ ' 
-,(l) r (l) r (2) 



(112a) 
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(112d) 



(112e) 



(112f) 



— i /ArWr^'r 1 

+ ^(8i(r 2 ;>) 3 ( r 2 2 >) 2 /i6) 
+ ^(i8rg ( r 2 :>) 2 ), 
= ^ 2 (63r«r 2 :> ( rg) 2 /2) 

+ / l 3 (243(rg) 3 (rg) 3 /32) 

+ /(63rg(r 2 ;>) 2 rg/2), 

= ^0 2 rf>) 

+ ^ 2 (54r« ( r 2 :>) 2 (r< 2 ») 2 ) 

+ a 3 (r 2 ^(729(r 2 |)) 3 (rg) 3 + 5i2r$rg)/64) 

+ ^(siigctS) 3 ^) 
+ ^(4(r 2 ;>) 2 (2ir< 2 » + 5r< 2 2 rg + 6or ( _ 2 » 4 

44I%)/5), 

= a 1 (6(r 2 2 >) 2 (2irgr 2 ;> + 6or« 4 rg + 35(r$) 

+4r -3,3)/ 5 ) 

+ ^ 3 (i053r^(rg) 2 (rg) 3 /8) 

+ /(3r 2 ^rg(729(r 2 ^) 3 (rg) 3 + 5632r^r®)/i28) 

+ ^ 5 (i053r®(r 2 ^) 3 (rg) 2 /8) 

+ A b (6(rg) 2 (2irgrg + 6or^ 4 rg + 35(rg) 2 

+4r (2 ] 3 )/5). (112h) 



The results are identical to those obtained by method of re- 
flections in Eqs. (lA27al >. (IA27bb . and dA27cl > for the terms 
containing one or two T's. The results reduce to those by Jef- 
frey and Onishi 13 in the no-slip limit y = and those by Keh 
and Chen 16 in the case of yi = y^. 



(H2g) 

2 



2. Y B Functions 

The problem for Y B is exactly the same as for Y A . The 
difference is that the force is calculated in Y A while the torque 
in Y B . Correspondingly, The coefficient fY Ea is defined as 



h 



2 2' 



2 ^Kk-qyq^' 



(113) 



3=0 



for Qi pt] obtained by the recurrence relations for Y A . The ex- 
plicit forms up to k = 7 are 
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cYB\ 



fl 



= 0, (114a) 

= 0, (114b) 

= A(-6T^rg), (114c) 

= iWW. (H4d) 

= ^ 2 (-27r^r^(rg) 2 /2), (114e) 

= ^(-i20 2 r 2 2 ))' 

+ ^ 2 (-8ir(;> ( r 2 ;») 2 (r 2 2 ») 2 /4) 

+ a 3 ^^^), ( ii4f) 

= ^(-losO^rg) 2 ) 

+ ^ 3 (-243r^ ( r^) 2 (rg) 3 /8) 

+ /(^rgrgrg), (U4 g) 

= ^ 2 (-i89(r[;j) 2 r^ ( rg) 2 ) 

+ a 3 (-3r^(256or[^r® + 243(r 2 |>) 3 (rg) 3 )/i6) 

+ ^(-243r^r< 2 ) ( r«) 2 rg) 

+ ^ 5 (48T^rg(7rf > - 6Tfl - 4r <2 1 ) 4 )) . ( 1 14h) 



The results are identical to those obtained by method of re- 
flections in Eqs. (IA32ab . (IA32bl ). and dA32cl ) for the terms 
containing one or two F's. The results reduce to those by Jef- 
frey and Onishi 13 in the no-slip limit y = and those by Keh 
and Chen 16 in the case of yj = 72 • 



3. Y Function 

With the same recurrence relations and the initial condition 
for Y A , that is, for the translating particles, the function Y G is 
obtained from the coefficient Pz pq for the stresslet instead of 
P\pq f°i" the force. 

In this case, the coefficient f7 Ga is defined as 



YCa 
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-12*' V P (a) A q 

■ q=0 



(115) 
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The explicit forms up to k — 7 are 
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(116h) 



The results are identical to those obtained by method of reflec- 
tions in Eqs. ( IA38at , ( IA38bb . and ( IA38cb for the terms con- 
taining one or two F's. The results reduce to those by Jeffrey 14 
in the no-slip limit y — 0. 



4. Y c Function 



The function Y c gives the torque for the rotating particles 
with m = 1. Therefore, it is derived by the coefficient Qi pq for 
the torque from the same recurrence relations as for Y A , but 
with different initial condition 
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In this case, the coefficient f? Ca is defined as 



f YCa _ 7 k y n (.a) ,q+j 
1=0 



(118) 



up to k — 7 are 
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(119h) 



The results are identical to those obtained by method of reflec- 
tions in Eqs. dA49a| ) and dA49bl > for the terms containing one 
or two T's. The results reduce to those by Jeffrey and Onishi 13 
in the no-slip limit y = and those by Keh and Cher. 16 in the 
case of ~y\ = j2. 



5. Y H Function 

With the same recurrence relations and the initial condition 
for Y c , that is, for the rotating particles, the function Y H is 
obtained from the coefficient p2 Pq for the stresslet instead of 
Qipq for the torque. 

In this case, the coefficient f™ a is defined as 
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' q=0 



) 

2(k-q)q 



(120) 



where j — for even k and j — 1 for odd k. The explicit forms 
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a 3 (-4or^rg(5rg - 21$)) , (i2i g ) 



up to k — 7 are 
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\ JVL OX 0.5 1 2.3 1 2,3/ 

+ ^(i80(r< 2 ') 2 r^r^). 



(121h) 



where j — for even k and j = 1 for odd k. The explicit forms 



The results are identical to those obtained by method of re- 
flections in Eqs. ( IA55al ) and dA55bt for the terms containing 
one or two T's. The results reduce to those by Jeffrey 14 in the 
no-slip limit y = 0. 



15 



Y M Function 



D. Z Functions (m = 2) 



The function Y M gives the stresslet under a strain flow for 
m = 1. Therefore, it is derived by the coefficient P 2pq for the 
stresslet from the same recurrence relations as for Y A , but with 
different initial condition. The boundary conditions are 



Xml = -^(-^) a a a E a 6i m 8 2n , 

mn 

Q) (a) = 



(122a) 



= ^-lfa a E a (l-6y)6 hn 6 2 n, (122b) 

(122c) 



which correspond to Eq. (J-54) with the correction due to the 
slip. The expansions used here are [in (J-55), (J-56), and (J- 
57)] 
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.,(<>') 



10 - V V YnEi J a 



*> = (-ly^-a, a E a > > nm \ "\X t\ , (123b) 
v > 3 "Z-jZ-j 2(2n+ 1) ff 3 -«' v ' 



(a) _\l^ a p \ \ n ff 

p=0 q=0 



10 



CO CO 



(123c) 



The boundary conditions are given by 
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<k = o, 



(127c) 



which correspond to Eq. (J-69) with the correction due to the 
slip. The expansions used here are [in (J-70), (J-71), and (J- 
72)] 
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From Eqs. ( I87ab . d87bb . and (I87cl > for m = 2 and the expan- 
sions above, the recurrence relations are given by 



The initial conditions are given from Eqs. (186a) , (186b) , and 
(l86cT > by 



> m -S«2l¥l vM=6 n2 T™ 0^ = 0. (124) 
In this case, the coefficient f™ a is defined as 

fI M " = ^f J P%-^ q+i , (125) 

q=Q 

where j — for even k and /' = 1 for odd k. The explicit forms 
up to k - 7 are 
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The results are identical to those obtained by method of reflec- 
tions in Eqs. ( IA67ab , (IA67bb . and (IA67ct for the terms con- 
taining one or two F's. The results reduce to those by Jeffrey 14 
in the no-slip limit y = 0. 



The initial conditions are obtained from Eqs. (I86ab . (I86bb , and 
(I86cb as 

P ( 2 = 62nT 2 % V«=WtS. G^o = 0. (130) 
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In this case, the coefficient ff Ma is defined as 
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The results are identical to those obtained by method of reflec- 
tions in Eqs. ( IA72ab . ( IA72bl ). and ( IA72cb for the terms con- 
taining one or two T's. The results reduce to those by Jeffrey 14 
in the no-slip limit y = 0. 



V. CONCLUDING REMARKS 

We have extended the calculations of resistance functions 
of two spheres with arbitrary size by the method of twin multi- 
pole expansions in general linear flows by Jeffrey and Onishi 13 
and Jeffrey 14 to the slip particles with the Navier slip boundary 
condition with arbitrary slip lengths. This extension comple- 
ments the previous results of slip particles obtained by Keh 
and Chen 16 for the same scaled slip lengths without strain 
flow. In limiting cases, the present calculations recover the ex- 
isting results, that is, those by Jeffrey et al.-^ill in the no-slip 
limit, and those by Keh and Chen 16 in the case of equal scaled 
slip lengths. We have also derived the resistance functions 
by the method of reflections and demonstrated its consistency 
with the twin multipole expansions. 

The present solutions of two-sphere problem cover much 
wider range than the previous solutions. Because the particle 
radii and slip lengths can be chosen independently, the solu- 
tions are not only applicable to the problem of two bubbles 
(demonstrated in Keh and Chen) 16 but also to that of solid 



particle and gas bubble, for example, with arbitrary sizes. In 
addition to these fundamental aspects in fluid dynamics, the 
solutions of slip particles is quite important for applications 
to micro- and nanofluidics, where the no-slip boundary con- 
dition may break. 14-6 Furthermore, the importance of the ex- 
act solution should be emphasized, because of the fact that 
the slip boundary condition is solved under relatively limited 
cases compared to the no-slip case. 

Using the multipole expansions and Faxen's laws derived in 
the present paper, recently the Stokesian dynamics method 26 
is extended from the no-slip particles to the slip particles. 12 
Because the lubrication corrections are missing in the formu- 
lation, the applicability is limited to relatively dilute configu- 
rations. The present work is a first step to improve the Stoke- 
sian dynamics method for slip particles at the level of the no- 
slip particles. To complete the program, we have to obtain the 
asymptotic forms of resistance functions by lubrication the- 
ory. To the authors' knowledge, just a few functions 10 are ob- 
tained for slip particles by now. On the other hand, the present 
exact solution expressed by 1 jr expansion is the complete set 
for the motion of rigid (slip) particles, that is, it contains all 
1 1 scalar functions for each pair of particles a/3, so that it is 
quite helpful to complete the lubrication theory for slip parti- 
cles and to develop the Stokesian dynamics method with lu- 
brication effect for arbitrary slip particles. 

The computer programs used in the paper and the results of 
coefficients for higher orders (up to k = 20) are available on 
the open source project "RYUON-twobody". 27 



Acknowledgments 

This work was supported by the National Research Council 
(NRC) of Canada. One of the authors (KI) thanks Professor 
David Jeffrey for his kind support and fruitful discussions. 



Appendix A: Method of Reflections 

Here we summarize the results of lower coefficients ob- 
tained by the method of reflections functions. 

1. Faxen's Laws 



From Eqs. (I50t . d56t . and d66l l in the previous section, 
the disturbance velocity field at position x caused by a sin- 
gle sphere a at x a with slip length y a is given by 

,2 



v(x) = 
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(i + r<^v 2 )j(x-xj-FW 
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p(a) 



1 + my a 
1 + ny a ' 



(Al) 



(A2) 



17 



and the force F^ a \ torque T <ff) , and stresslet S^ a) on the sphere 
are given by 
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(A3a) 
(A3b) 

(A3c) 



(See Eqs. d49i >, (l55l >. and (l63l in the previous section.) Read- 
ing Eq. ( IA1I ) as multipole expansion of the velocity field, 
Faxen's laws for slip sphere are derived as 
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where u' is the velocity field in absent of particle a. For later 
use, we rewrite Eq. (IAU in the resistance form by replacing 
F (a \ T ia \ and S (a) by U {a \ fl {a \ and from Eqs. (TAlal l. 
jAlbl , and dA3c] i as 
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if (a; - x„) : E (a \A7) 



2. Translating Spheres in Axisymmetric Motion 

Here we set the relative vector between particle 1 and 2 in 
Z direction as 



r = X2 - x\ — (0, 0, r). 



(A8) 



For the function X A , we set the velocity of the particle 1 par- 
allel to r as 



i7 (1) =(0,0, t/ (1) ). 



(A9) 



From Faxen's law for the force (IA41 with the disturbance 
field (IA7b with Eq. 1A91 , we have the force on the particle 2 
due to the translating particle 1 as 
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In terms of the scalar functions X A „, the force is expressed as 

Ff } = 6n[ia 2 X A 2 (s,A)U (2) S iz 

+37ryu(a 2 + ai)X^(5,/l)C/ (1) ^, (All) 
where s and A are defined in Eq. ®. Therefore, 
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From the symmetry of X A ^ in Eq. ITal i. we have 
-2 
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From the expression of X A 2 in Eq. ( I68bb . we have f* A as 



fXA _ o r (l) r (2) -, 

Jl ~ J1 2,3 1 2,3 A ' 

fXA _ Ajtr^r^ _zu 3 r (1 )r( 2) 

J3 ~ 0,3*2,3 H/tl 2,3 1 0,3- 



For the self part X A { , we have 
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These coefficients (and those for the rest of the functions be- 
low) will be compared with the results by twin multipole ex- 
pansions in Sec. |IV A| 

From Faxen's law for the torque (1A5K we have torque on 
the particle 2 due to the translating particle 1 as 



r(2) 



= o, 



(A15) 



because Sl a — in the present problem and djUi is symmet- 
ric about the indices j, k. This fact reflects that there is no X B 
function in Eq. 16b) . 

From Faxen's law for the stresslet 
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In terms of the scalar functions X^, the stresslet is expressed 
as 

S{f = ^(fl 2 + fli) 2 ^^ 1 ' I^j, - Itfyj , (A17) 



so that 
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15 



(1 +A) 2 s 



_p(2) r (l) + 



60 



,2 r 2 2.5 1 2,3 



p(2) r (l) 



(1 + i)V 2 ' 5 °< 3 



36^ 2 



r (2) r (D 



(A18) 
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From the symmetry of in Eq. d7fl ). we have 
-4 



vC _ 

12 " (l+A) 2 



(i + /i) 2 5 2 



(1 + /l)V 



From the expression of X^ 2 in Eq. ( I73bb . we have ff c as 



(A19) 



XG 




fi 



xg 



(A20a) 
(A20b) 



- 60/i 3 r^rg. (A20c) 



3. Translating Spheres in Asymmetric Motion 

Next, we study the asymmetric motion of the spheres to 
their center-to-center vector, that is, the velocity U m is in y- 
direction as 



J7 (1) = (0, t/ (1) ,0). 



Note that, for r = (0, 0, r), from Eq. d6at . we have 



1 aB u * 



1 a a u y 



(A21) 



(A22) 



From Faxen's law for the force (IA41 with the disturbance 
field dA7| ) with Eq. ( IA21I ). we have the force on the particle 2 
due to the translating particle 1 as 



U m 5 iy . (A23) 



1 flt 



_2 

- J2p(2) r (l) , ^2^]_ r (2) r (l) 
T 4 r 3 2,3 0,3 4 f 3 0,3*2,3 

In terms of the scalar functions Y A „, the force is expressed as 



7 (2) 



6^ua 2 Ff 2 (i,>l)C/ (2) ^ 

+37TjU(fl 2 + fll)l2l( 5 ' 



(A24) 



Therefore, 



y*(M) = rg 



2i (3 



1 



<2) r (D 



1 +/l\2(l +i)s 



F ' r 



2,3 X 2,3 



■(iSlg+^Slg) • (A25b) 



(1 +^) 3 s 3 

From the symmetry of Y A „ in Eq. ( TTbt , we have 



J?i(M) = r«, 



(A26a) 



Ff 2 (*,A) = - 



r (l)p(2) 

1 +/l\2(l + 2 ' 3 2 > 3 



(1 +/l)V 



(A 3 r«rg + <rg) . ( A26b) 



From the expression of Y A 2 in Eq. d69b| i, we have f7 A as 

fl A = 4« (A27a) 



ft A _ ^ r (l) r (2) , 
J I - -) 1 2.3 1 2.3 /l ' 



ft 



2 2,3* 2,3' 

iyWU 3 + iyW\a 



2,3 0,3 



0.3 1 2,3 



(A27b) 
(A27c) 



From Faxen's law for the torque (lA5t . we have the torque 
on the particle 2 due to the translating particle 1 as 



T (2) = -6n^a^V^T 2 ^U w 5 ix . (A28) 
In terms of the scalar functions YjL, the torque is expressed as 
T {2) = 4n^alY% 2 Si X U i2) + Mai +a l ) 2 Y 2 3 l 6 ix U (l) (A29) 



Therefore, 



Y22 = 0, 



Y B - 

I 2\ ~ 



-4 



6A 3 



<2) r (l) 



r'-'r 



(l+A) 2 (l+A) 2 s 2 °- 3 2 - 3 ' 
From the symmetry of Y B „ in Eq. dTcj i, we have 



Y B - 
1 12 _ 



-6/1 



(1 +/1) 2 (1 +/l) 2 S 



r (l) r (2) 



2 r 2 0,3' 2,3 



(A30a) 
(A30b) 

(A3 la) 
(A3 lb) 



From the expression of Y B { and Y B 2 in Eqs. ( |70at and d70b| 



we have fJ B as 



/r = 0, 

rYB 
J2 



_ 6/}r (l) r (2) 
U/11 0.3 1 2,3- 



From Faxen's law for the stresslet 
given by 



(A32a) 
(A32b) 
(A32c) 

the stresslet is 



S (2) = 



-f^ a 2 r 2 2 5 



'3fl 3 

A? 
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Lr (1) 1 r (2) fl 2 9fl i r (i)| 

4r 4 1 0,3 +1 0,2 102r 4 i 2,3 



(A25a) Note that 



S {2) = 



Anixa\Gf jk vf + nfi(a 2 + a^G^U^, 



(A33) 



(A34) 



where 



for e = (0, 0, 1) and U = (0, U, 0). Therefore, we have 



(A35) 



'22 

yG 
r 21 



0. 



4/1 



.(2) r (l) 



20A 2 

(l+A) 2 \(1 +A) 4 s^ 2 ^ ' 3 
12 A 3 



5 (1 +A) 4 s- 



_ r (2) r (l) 

! 4 C 4 1 0,5 1 2,3 



(A36a) 



(A36b) 
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From the symmetry of Y G in Eq. ( 7g I, we have 



4 / 2(M 3 



12 



p(l) r (2) 

(l+A) 2 \(l+A) 4 s 4 2 - 5 °- 3 



b. Torque in Asymmetric Motion For the asymmetric 
motion to the center-to-center vector, from Faxen's law for the 
(A37a) torque dASb with the disturbance field dATt with Eq. (IA39bb . 

we have the torque on the particle 2 due to the translating par- 
ticle 1 as 



12 A 



r (l) r (2) 



(1+A)V °- 5 2 - 3 
From the expression of Y^ 2 in Eq. ( I74bb . we have fJ G as 



(A38a) 

- w. (A38b) 

// G = 20/i 3 r^r^ + i2/ir^r^. (A38c) 



4. Rotating Spheres 

Next, we consider rotating spheres. In the two-body prob- 
lem with r = (0, 0, r), we set the angular velocity f2 (1) for the 
axisymmetric case by 

nf } = n (V) 6 iz , (A39a) 
and for the asymmetric case to the axis r by 



(A39b) 



a. Torque in Axisymmetric Motion From Faxen's law 
for the torque (IA5) with the disturbance field (IA7t with Eq. 
( IA39ab . we have the torque on the particle 2 due to the trans- 
lating particle 1 as 



(A37b) if = STr^rgfif + Anpft™ 1$ ^<^ (1) ■ (A45) 

In terms of the scalar functions F£j, the torque is expressed as 
Tf ] = 87r/ia^F 2 c 2 fi (2) (5, v +n/j(a 2 + atfY 2 = l Q. m S iy (A46) 



Therefore, 



yC _ r (2) 
2 22 _ i 0,3' 



Y c - 

I 2\ ~ 



4^ 3 



_ r (2) r 0). 

. A T A A 



(i + ^) 3 °> 3 °- 3 (i + w 

From the symmetry of YfL in Eq. ( TTel l. we have 



(A47a) 
(A47b) 



<l) r (2) 



_r u 'r 
r a a o a i 



8i 3 



(1+/1) 3 i 0,3 i 0,3 (1+/i) 3 i 3- < A48 > 



From the expression of Yf. in Eq. ( I72bb . we have /, FC as 



/r - o, 

f YC _ /ii3 r (l) r (2) 
J 3 0.3 1 0,3- 



c. Stresslet Because <9/z4 for the axisymmetric motion 
is anti-symmetric for z and j, there is no contribution to the 
stresslet. For the asymmetric motion to the axis, from Faxen's 
law for the stresslet (IA6t with the disturbance field (1A7K 



(A49a) 
(A49b) 



Tf ] = 87rjt/a 3 r^Q| 2) - Zn^alrV^T ( ^5 iz £l m . (A40) 



In terms of the scalar functions X^, the torque is expressed as 
Tf> = 2,TiiialX^ 2 Q. {1) Si z + nix{a2 + a l fX^Q. m 5 iz (Ml) 



Therefore, 



yC _ r (2) 
A 22 ~ 1 0,3' 



(A42a) 



A 21 - 



8A 3 8 



rgr$. (A42b) 



(1 + /T) 3 (1 +^) 3 5 3 

From the symmetry of in Eq. d7dl >. we have 



8 8,i 3 



XhfA) = — r — " , ^ hfl (A43) 

12 (1 + T) 3 (1 +A) 3 ,s 3 °- 3 °- 3 v ; 

From the expression of Xf 2 in Eq. d71bb . we have f* c as 



fXC _ r (i) 

■'0 _ 1 0,3' 

ff C = o, 

fXC _ oi3 r (l) r (2) 



(A44a) 
(A44b) 
(A44c) 



s ® = lQa^^igrJg (cv,() ;i + «)..,«),.) (A50) 

The stresslet on particle 2 caused by particle 1 is given by 



Mai + ayfH^Q.^, 
where, for r = (0, 0, r) and dp = Q. <l) 6 ky , 



Therefore, 



Y" - 



10 8 



,(2) r (l) 



(l+Af (1 +T) 3 s 3 2 - 5 °- 3 ' 
From the symmetry of Y^„ in Eq. dThl i, we have 



12 ~ (1 + ^) 3 (1 + i) 3 * 3 2 ' 5 °' 3 ' 



(A51) 



rf$nf = Y% (<M n + <V, ()/1 ) (A52) 



(A53) 



(A54) 



From the expression of YK in Eq. ( I75bl i. we have f7 H as 



fl H = o, 



(A55a) 



fl H = lO^r®. (A55b) 
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5. Spheres in Strain Flow 

Next, we consider the problem under the strain flow. Let us 
define three types of strain by 



El = E x 



(**-*). 



E Y U = E Y {S kz Six + $kx6h) , 
E z = E z (d kx 6 lx -d ky 6 ly ), 



(A56a) 

(A56b) 
(A56c) 



which correspond to the scalar functions X^, Y^, and Z^, 
respectively. 

In the following, we will see 5' (2;1) , the stresslet on parti- 
cle 2 caused by particle 1, which is related to the resistance 
functions X^ v Y^ 2 , and Z^ 2 . From Faxen's law for the stresslet 



it is given by 



c (2;l) 20 3 (2) 

y = y^ fl 2 r 2,5 



, r (2) "2 V 

°- 2 10 , 



Udiuf + djlip] (B 2 ) 



(A57) 

d. Function X M Substituting the disturbance field (lA7t 
with Eu ( IA56ab into Eq. (IA57b . we have 



<(2;1) 



20 3 

= -J nfxa l 



5a 3 



_J. r (2) r (l) 
„3 1 2,5*2,5 



— ( a 5 r {2) r (1) +aVr (2) r (1 M 

r 5 \ a l L 2,5 L 0,5 + a 2 a l l 0,5* 2,5) 
(5, 



x{6iz 6 jz -f)EX 



In terms of the scalar function XM, it is written as 



SftV = -ntfto + atfX™ \5 k 6 k - -f\E 



3 V M 



(A58) 



(A59) 



Therefore, 



X 
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(l+A) 3 

192 



40 



r (2) r (l) 



a + i)V 2 - 5 2 - 5 

/ r (2) r (l) i2 r (2) r (l)\ 
I 1 2.5 1 0,5 + A 1 0.5 1 2,5 j 



(1 +A) 5 S 5 

From the symmetry of X^ in Eq. ((TI]), we have 

40^ 3 



(A60) 
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(1+^) 3 
192 



r (l) r (2) 

i M CM 



(l+i)V 2 ' 5 2 ' 5 

I ,5 r (l) r (2) 3 r (l) r (2)\ 
1 2.5 1 0,5 + A 1 0.5 1 2,5 j 



(1 +T)V 

From the expression of X^ 2 in Eq. ( I76bb . we have /* M as 



(A61) 



■XM 



(A62a) 
(A62b) 

I92(^rg + ^ r «rg). ( A62c) 



40^F«Fg 



e. Function Y M Substituting the disturbance field ( IA7b 
with Ej, ( IA56bb into Eq. ( IA57b . we have 



.-.(2;!) 



20 3 
= —njia 2 



_5 «i_ (2) (1) 
2 r 3 1 2,5 1 2,5 



h lf fl 5 r (2) r (D 2 3 r (2) r (l)\ 
h r 5 \ a l L 2,5 L 0,5 + fl 2 fl l 1 0.5 1 2,5 J 

i(d ix 6 jz + d b 8 jx )E Y . 



(A63) 



In terms of the scalar function Y¥., it is written as 



Sfr 1} = -nn{a 2 + a^Y™ {d iz 6 jx + 6 ix 6 jz ) E Y . (A64) 



Therefore, 
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128 
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(1 + A) 5 s 5 

From the symmetry of Y% in Eq. (|7j]l, we have 

20/1 3 



(A65) 
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/ ,5 r (l) r (2) j3 r (l) r (2)\ 
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(1 + A) 5 S 5 

From the expression of Y^ 2 in Eq. d77bb . we have f YM as 



(A66) 



/r = 0, 



/ 3 ™ = -20^T^T 25 . 



3 r (l) r (2) 



(A67a) 
(A67b) 

/™ = 128 (^rgrg + ^ 3 r^rg) . (A67c) 

/ Function Z M Substituting the disturbance field ( IA7b 
with E z ( IA56cl ) into Eq. ( IA57I ). we have 



:(2;i) 
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7T//fl 2 ^5 ^l 1 2.5 1 0,5 + a 2 a l l 0.5 1 2,5 ) 

x(s ix S jx -S iy 8jy)E z . (A68) 
In terms of the scalar function Z 2{ , it is written as 

Sff ) = -n(x(a 2 + aifZ^ (6 ix 6 jx - 5 iy 6 h ) E z . (A69) 



Therefore, 



Z* = -I 



32 



. -(r^ + ^rWn. (A70) 

(1 +i) 3 (1 +A) 5 s 5 \ 2 - 5 °- 5 °- 5 2 ' 5 > 
From the symmetry of Z^„ in Eq. d7ki . we have 
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^12 _ /i , ,-,3 
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From the expression of Z^ in Eq. ( I78bb . we have f ZM as 



fi 



■ZM 



(A72a) 
(A72b) 



f ZM = 32(^r^r[ 2 > + i 3 r^rg). (A72c) 
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